# Load Packages
library(XLConnect)
library(xtable)
library(scales)
library(stargazer)

av_interaction_flag <- TRUE
choice_seq_flag <- FALSE
all_rs_variables <- FALSE

# Make Zero Residual Counterfactual Tables			

	setwd("C:/Users/eas24f/OneDrive - Florida State University/Research/Inertia RESTAT/Data")
	
	oligopoly_scenarios <- c("ACA","Base","No_Inertia_NoSupply","No_Inertia","No_RA","No_Inertia_No_RA")
	pc_scenarios <- c("PC","No_Inertia_PC","No_RA_PC","No_Inertia_NO_RA_PC")
	if(all_rs_variables) {
		pc_scenarios <- c("PC","No_Inertia_PC")
	}
	norapc_scenarios <- c("No_RA_PC","No_Inertia_NO_RA_PC")
	inertia_columns <- c(oligopoly_scenarios,pc_scenarios)
	year_columns <- c("Base 2015","No Inertia 2015","Base 2016","No Inertia 2016","Base 2017","No Inertia 2017","Base 2018","No Inertia 2018")
	
	example_data <- read.csv("2015Base0_supply_output_intvars_avint.csv",header=TRUE,row.names=1)
	fields <- rownames(example_data)
	cw_fields <- c("CWlt250","CW250to400","CWgt400","CW0to17","CW18to34","CW35to54","CW55plus",
						"CWFemale","CWMale","CWSingle","CWFamily","CWAsian","CWBlack","CWHispanic","CWWhite","CWOther")
	other_fields <- setdiff(fields,cw_fields)
	data_inertia <- as.data.frame(matrix(NA,length(fields),length(inertia_columns),dimnames=list(fields,inertia_columns)))
	data_inertia_years <- as.data.frame(matrix(NA,length(other_fields),length(year_columns),dimnames=list(other_fields,year_columns)))
	
	for(scenario in inertia_columns) {
		if(all_rs_variables) {
			data2015 <- read.csv(paste("2015",scenario,"0_supply_output_allvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			data_inertia[other_fields,scenario] <- data2015
			if(scenario == "ACA") {
				data2015base <- read.csv(paste("2015",scenario,"0_supply_output_intvars_seq_avint.csv",sep=""),header=TRUE,row.names=1)
			}
		} else if(choice_seq_flag & av_interaction_flag) {
			data2015 <- read.csv(paste("2015",scenario,"0_supply_output_intvars_seq_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2016 <- read.csv(paste("2016",scenario,"0_supply_output_intvars_seq_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2017 <- read.csv(paste("2017",scenario,"0_supply_output_intvars_seq_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2018 <- read.csv(paste("2018",scenario,"0_supply_output_intvars_seq_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			data_inertia[,scenario] <- (data2015 + data2016 + data2017 + data2018)/4
			if(scenario %in% norapc_scenarios) {
				data_inertia[,scenario] <- (data2015 + data2016 + data2017)/3
			}
			data_inertia[,scenario] <- data2015
		} else if (av_interaction_flag) {
			data2015 <- read.csv(paste("2015",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			
			if(scenario %in% c("Base","No_Inertia")) {
				data_inertia[,scenario] <- data2015
			} else {
				data_inertia[other_fields,scenario] <- data2015
			}
			
			if(scenario == "ACA") {
				data2015base <- read.csv(paste("2015",scenario,"0_supply_output_intvars_seq_avint.csv",sep=""),header=TRUE,row.names=1)
			} else if(scenario == "Base") {
				data_inertia_years[,"Base 2015"] <- read.csv(paste("2015",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[other_fields,1]
				data_inertia_years[,"Base 2016"] <- read.csv(paste("2016",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
				data_inertia_years[,"Base 2017"] <- read.csv(paste("2017",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
				data_inertia_years[,"Base 2018"] <- read.csv(paste("2018",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			} else if(scenario == "No_Inertia") {
				data_inertia_years[,"No Inertia 2015"] <- read.csv(paste("2015",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[other_fields,1]
				data_inertia_years[,"No Inertia 2016"] <- read.csv(paste("2016",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
				data_inertia_years[,"No Inertia 2017"] <- read.csv(paste("2017",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
				data_inertia_years[,"No Inertia 2018"] <- read.csv(paste("2018",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			}
		} else {
			data2015 <- read.csv(paste("2015",scenario,"0_supply_output_noseq.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2016 <- read.csv(paste("2016",scenario,"0_supply_output_noseq.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2017 <- read.csv(paste("2017",scenario,"0_supply_output_noseq.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2018 <- read.csv(paste("2018",scenario,"0_supply_output_noseq.csv",sep=""),header=TRUE,row.names=1)[,1]
			data_inertia[,scenario] <- (data2015 + data2016 + data2017 + data2018)/4
		}
	}
	
	
	if(av_interaction_flag & !choice_seq_flag & !all_rs_variables) {
		data_inertia["Profit",setdiff(inertia_columns,"ACA")] <- data_inertia["Profit",setdiff(inertia_columns,"ACA")] -
			(data_inertia["Profit","ACA"] - data2015base["Profit",1])
		data_inertia_years["Profit",] <- data_inertia_years["Profit",] - (data_inertia["Profit","ACA"] - data2015base["Profit",1])
	}
	data_inertia["Profit",pc_scenarios] <- 0
		
	# ACA/Base Metal Coverage
	if(av_interaction_flag) {
		nonsilver_plans <- c("Bronze","Gold","Platinum")
		silver_error <- data_inertia["Silver","ACA"] - data_inertia["Silver","Base"]
		
		data_inertia["Silver",setdiff(inertia_columns,"ACA")] <- data_inertia["Silver",setdiff(inertia_columns,"ACA")] + silver_error
		data_inertia[nonsilver_plans,setdiff(inertia_columns,"ACA")] <- 
			data_inertia[nonsilver_plans,setdiff(inertia_columns,"ACA")] - t(t(data_inertia[nonsilver_plans,setdiff(inertia_columns,"ACA")])/
			colSums(data_inertia[nonsilver_plans,setdiff(inertia_columns,"ACA")])) * silver_error
			
		#data_inertia_years["Silver",] <- data_inertia_years["Silver",] + silver_error
		#data_inertia_years[nonsilver_plans,] <- data_inertia_years[nonsilver_plans,] - 
		#	t(t(data_inertia_years[nonsilver_plans,])/colSums(data_inertia_years[nonsilver_plans,])) * silver_error
		
	}
	
	# Apply deadweight loss of taxation (from DPR AER Paper) to government spending
	data_inertia["Social Welfare",] <- data_inertia["Consumer Surplus",] + data_inertia["Profit",] - 
		1.3 * (data_inertia["Premium Subsidies",] + data_inertia["Cost Sharing Subsidies",] - 
			data_inertia["Penalties",] + data_inertia["Uncompensated Care",])
	data_inertia_years["Social Welfare",] <- data_inertia_years["Consumer Surplus",] + data_inertia_years["Profit",] - 
		1.3 * (data_inertia_years["Premium Subsidies",] + data_inertia_years["Cost Sharing Subsidies",] - 
			data_inertia_years["Penalties",] + data_inertia_years["Uncompensated Care",])
	
	# Year adjustment_factor
	data_inertia_years[,c("Base 2015","Base 2016","Base 2017","Base 2018")] <- data_inertia_years[,c("Base 2015","Base 2016","Base 2017","Base 2018")] * 
		data_inertia_years[,"Base 2015"] / rowMeans(data_inertia_years[,c("Base 2015","Base 2016","Base 2017","Base 2018")])
	data_inertia_years[,c("No Inertia 2015","No Inertia 2016","No Inertia 2017","No Inertia 2018")] <- data_inertia_years[,c("No Inertia 2015","No Inertia 2016","No Inertia 2017","No Inertia 2018")] * 
		data_inertia_years[,"No Inertia 2015"] / rowMeans(data_inertia_years[,c("No Inertia 2015","No Inertia 2016","No Inertia 2017","No Inertia 2018")])
	
	
	if (choice_seq_flag) {
		write.csv(data_inertia,"averaged_supply_output_inertia_zeroresid_intvars_seq_avint.csv")
	} else if (all_rs_variables) {
		write.csv(data_inertia,"averaged_supply_output_inertia_zeroresid_allvars_avint.csv")
	} else if(av_interaction_flag) {
		write.csv(data_inertia,"averaged_supply_output_inertia_zeroresid_intvars_avint.csv")
		write.csv(data_inertia_years,"years_supply_output_inertia_zeroresid_intvars_avint.csv")
	} else {
		write.csv(data_inertia,"averaged_supply_output_inertia_zeroresid_noseq.csv")
	}
	
	# Extract Output
	prem_rows <- c("Bronze Premium","Silver Premium","Gold Premium","Platinum Premium",
		"Anthem Premium","Blue Shield Premium","Health Net Premium","Kaiser Premium","Other Insurer Premium",
		"HMO Premium","PPO Premium","Average Premium","Average Subsidized Premium")
	coverage_rows <- c("Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"% Switching Plan (Base)","% Switching Metal (Base)","% Switching Insurer (Base)",
		#"% Switching Plan (Setting)","% Switching Plan (Setting)","% Switching Plan (Setting)"
		"Anthem","Blue_Shield","Health_Net","Kaiser","Small_Insurer")
	#coverage_rows <- c("Bronze","Silver","Gold","Platinum","Total Coverage","Percent Enrolled","Percent Switching")
	coverage_names <- c("Bronze","Silver","Gold","Platinum")
	coverage_names_insurers <- c("Anthem","Blue Shield","Health Net","Kaiser","Other Insurer")
	claims_rows <- c("Bronze Avg Claims Resid","Silver Avg Claims Resid","Gold Avg Claims Resid","Platinum Avg Claims Resid","Average Claims Resid")
	welfare_rows <- c("Consumer Surplus","Profit","Premium Subsidies","Cost Sharing Subsidies",
				"Penalties","Uncompensated Care","Social Welfare")
	welfare_names <- c("Cons. Surplus","Profit","Prem. Subsidies","CSRs",
		"Penalties","Uncomp. Care","Social Welfare")
				
	output <- data_inertia[c(prem_rows,coverage_rows,claims_rows),setdiff(inertia_columns,"ACA")]
	rownames(output) <- c("Bronze P","Silver P","Gold P","Platinum P",
		"Anthem P","Blue Shield P","Health Net P","Kaiser P","Other Insurer P","HMO","PPO","Average","Average (w/ Subsidies)",
		"Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"Plan","Metal","Insurer",
		"Anthem","Blue Shield","Health Net","Kaiser","Other Insurer",
		"Bronze Cl","Silver Cl","Gold Cl","Platinum Cl","Average Cl")
	colnames(output) <- c("Base","(1)","(2)","(3)","(4)","(5)","(6)","(7)","(8)")
	
	report_shares <- TRUE
	if(report_shares) {
		for(s in 1:dim(output)[2]) {
			output[setdiff(coverage_names,"Total Coverage"),s] <- 100* output[setdiff(coverage_names,"Total Coverage"),s]/output["Total Coverage",s]
			output[setdiff(coverage_names_insurers,"Total Coverage"),s] <- 100* output[setdiff(coverage_names_insurers,"Total Coverage"),s]/output["Total Coverage",s]
		}
	}
	
	output[c("% Enrolled","Plan","Metal","Insurer"),] <- output[c("% Enrolled","Plan","Metal","Insurer"),] * 100
	
	indent_tags <- c(FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE)
	rownames(output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(output)[indent_tags],sep="")
	
	print(xtable(output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
	# Create welfare table separately
	
	welfare_output <- data_inertia[welfare_rows,setdiff(inertia_columns,"ACA")] * 12
	rownames(welfare_output) <- welfare_names
	colnames(welfare_output) <- c("Base","(1)","(2)","(3)","(4)","(5)","(6)","(7)","(8)")
	#welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1]) * 12
	welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1])
	indent_tags <- c(FALSE,FALSE,TRUE,TRUE,TRUE,TRUE,FALSE)
	rownames(welfare_output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(welfare_output)[indent_tags],sep="")
	
	print(xtable(welfare_output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
	# Appendix Welfare Table	

	welfare_rows <- c("Consumer Surplus",
				"Consumer Surplus - Kappa 75%","Consumer Surplus - Kappa 50%","Consumer Surplus - Kappa 25%","Consumer Surplus - Kappa 0%",
				"Social Welfare",
				"Social Welfare - Kappa 75%","Social Welfare - Kappa 50%","Social Welfare - Kappa 25%","Social Welfare - Kappa 0%")
	welfare_names <- c("Cons. Surplus","Kappa 75%","Kappa 50%","Kappa 25%","Kappa 0%",
		"Social Welfare","SW Kappa 75%","SW Kappa 50%","SW Kappa 25%","SW Kappa 0%")
		
	output <- data_inertia[welfare_rows,setdiff(inertia_columns,"ACA")]
	rownames(output) <- welfare_names
	colnames(output) <- c("Base","(1)","(2)","(3)","(4)","(5)","(6)","(7)","(8)")
	
	#output[welfare_names,] <- (output[welfare_names,] - output[welfare_names,1]) * 12
	output[welfare_names,] <- (output[welfare_names,] - output[welfare_names,1]) 
	output["SW Kappa 75%",] <- output["Social Welfare",] + (output["Kappa 75%",] - output["Cons. Surplus",])
	output["SW Kappa 50%",] <- output["Social Welfare",] + (output["Kappa 50%",] - output["Cons. Surplus",])
	output["SW Kappa 25%",] <- output["Social Welfare",] + (output["Kappa 25%",] - output["Cons. Surplus",])
	output["SW Kappa 0%",] <- output["Social Welfare",] + (output["Kappa 0%",] - output["Cons. Surplus",])
		
	indent_tags <- c(FALSE,TRUE,TRUE,TRUE,TRUE,
					FALSE,TRUE,TRUE,TRUE,TRUE)
	rownames(output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(output)[indent_tags],sep="")
	
	print(xtable(12*output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
	# Create Premium Selection Table

		# Extract output
		prem_rows <- c("Bronze Premium (No Selection)","Silver Premium (No Selection)",
							"Gold Premium (No Selection)","Platinum Premium (No Selection)","Average Premium (No Selection)",
							"Bronze Premium (No Metal Selection)","Silver Premium (No Metal Selection)",
							"Gold Premium (No Metal Selection)","Platinum Premium (No Metal Selection)","Average Premium (No Metal Selection)",
							"Bronze Premium (Intensive Selection)","Silver Premium (Intensive Selection)",
							"Gold Premium (Intensive Selection)","Platinum Premium (Intensive Selection)","Average Premium (Intensive Selection)",
							"Bronze Premium","Silver Premium","Gold Premium","Platinum Premium","Average Premium")
		
		prem_output <- data_inertia[prem_rows,setdiff(inertia_columns,"ACA")]
		rownames(prem_output) <- c("Bronze P","Silver P","Gold P","Platinum P","Average P",
			"Bronze P2","Silver P2","Gold P2","Platinum P2","Average P2",
			"Bronze P3","Silver P3","Gold P3","Platinum P3","Average P3",
			"Bronze P4","Silver P4","Gold P4","Platinum P4","Average P4")
		colnames(prem_output) <- c("Base","(1)","(2)","(3)","(4)","(5)","(6)","(7)","(8)")
		
		print(xtable(prem_output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
			type="latex",sanitize.text.function=function(x){x})	
	
	# Create Claims Selection Table

		# Extract output
		claims_rows <- c("Bronze Avg Claims (No Selection)","Silver Avg Claims (No Selection)",
							"Gold Avg Claims (No Selection)","Platinum Avg Claims (No Selection)","Average Claims (No Selection)",
							"Bronze Avg Claims (No Metal Selection)","Silver Avg Claims (No Metal Selection)",
							"Gold Avg Claims (No Metal Selection)","Platinum Avg Claims (No Metal Selection)","Average Claims (No Metal Selection)",
							"Bronze Avg Claims (Intensive Selection)","Silver Avg Claims (Intensive Selection)",
							"Gold Avg Claims (Intensive Selection)","Platinum Avg Claims (Intensive Selection)","Average Claims (Intensive Selection)",
							"Bronze Avg Claims Resid","Silver Avg Claims Resid","Gold Avg Claims Resid","Platinum Avg Claims Resid","Average Claims Resid")
		
		claims_output <- data_inertia[claims_rows,setdiff(inertia_columns,"ACA")]
		
		# Fix Base Column
		claims_output[,"Base"] <- claims_output[c("Bronze Avg Claims Resid","Silver Avg Claims Resid","Gold Avg Claims Resid","Platinum Avg Claims Resid","Average Claims Resid"),"Base"]
		
		# Fix No Selection Row
		claims_output[c("Bronze Avg Claims (No Selection)","Silver Avg Claims (No Selection)","Gold Avg Claims (No Selection)","Platinum Avg Claims (No Selection)","Average Claims (No Selection)"),] <- 
			claims_output[c("Bronze Avg Claims (No Selection)","Silver Avg Claims (No Selection)","Gold Avg Claims (No Selection)","Platinum Avg Claims (No Selection)","Average Claims (No Selection)"),"Base"]
		
		
		rownames(claims_output) <- c("Bronze C","Silver C","Gold C","Platinum C","Average C",
			"Bronze C2","Silver C2","Gold C2","Platinum C2","Average C2",
			"Bronze C3","Silver C3","Gold C3","Platinum C3","Average C3",
			"Bronze C4","Silver C4","Gold C4","Platinum C4","Average C4")
		colnames(claims_output) <- c("Base","(1)","(2)","(3)","(4)","(5)","(6)","(7)","(8)")
		
		print(xtable(claims_output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
			type="latex",sanitize.text.function=function(x){x})	
			
# Make Standard Error Tables
	# Two scenarios: Base and No Inertia
	# Fields: mean, std error, 2.5%, 97.5% (4 * 2 = 8 total)
			
	ns <- 100		
			
	output_columns <- c("Mean","St. Err","2.5%","97.5%","Mean NI","St. Err NI","2.5% NI","97.5% NI")		
	prem_rows <- c("Bronze Premium","Silver Premium","Gold Premium","Platinum Premium",
		"Anthem Premium","Blue Shield Premium","Health Net Premium","Kaiser Premium","Other Insurer Premium",
		"HMO Premium","PPO Premium","Average Premium","Average Subsidized Premium")
	coverage_rows <- c("Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"% Switching Plan (Base)","% Switching Metal (Base)","% Switching Insurer (Base)",
		#"% Switching Plan (Setting)","% Switching Plan (Setting)","% Switching Plan (Setting)"
		"Anthem","Blue_Shield","Health_Net","Kaiser","Small_Insurer")
	claims_rows <- c("Bronze Avg Claims Resid","Silver Avg Claims Resid","Gold Avg Claims Resid","Platinum Avg Claims Resid","Average Claims Resid")
	welfare_rows <- c("Consumer Surplus","Profit","Premium Subsidies","Cost Sharing Subsidies",
			"Penalties","Uncompensated Care","Social Welfare")
	welfare_names <- c("Cons. Surplus","Profit","Prem. Subsidies","CSRs",
		"Penalties","Uncomp. Care","Social Welfare")
				
	
	all_rownames <- c(prem_rows,coverage_rows,claims_rows,welfare_rows)
	output <- matrix(NA, length(all_rownames),length(output_columns),dimnames=list(all_rownames,output_columns))
	
	output[,"Mean"] <- data_inertia[all_rownames,"Base"]
	output[,"Mean NI"] <- data_inertia[all_rownames,"No_Inertia"]
	base_object <- noinertia_object <- matrix(NA,nrow(data_inertia),ns)
	rownames(base_object) <- rownames(data_inertia)
	rownames(noinertia_object) <- rownames(data_inertia)
	
	rownames(base_object)[which(rownames(base_object) == "Average Claims Resids")] <- "Average Claims Resid"
	rownames(noinertia_object)[which(rownames(noinertia_object) == "Average Claims Resids")] <- "Average Claims Resid"
	
	for(s in 1:ns) {
		base_object[,s] <- read.csv(paste("2015Base0_supply_output_intvars_avint_",s,".csv",sep=""),header=TRUE,row.names=1)[,1]
		noinertia_object[,s] <- read.csv(paste("2015No_Inertia0_supply_output_intvars_avint_",s,".csv",sep=""),header=TRUE,row.names=1)[,1]
	}
	
	for(n in all_rownames) {
		output[n,"St. Err"] <- sd(base_object[n,])/sqrt(ns)
		output[n,"2.5%"] <- quantile(base_object[n,],0.025)
		output[n,"97.5%"] <- quantile(base_object[n,],0.975)
		output[n,"St. Err NI"] <- sd(noinertia_object[n,])/sqrt(ns)
		output[n,"2.5% NI"] <- quantile(noinertia_object[n,],0.025)
		output[n,"97.5% NI"] <- quantile(noinertia_object[n,],0.975)
	}
	
	
	output["Profit",c("2.5%","97.5%")] <- c(output["Profit","Mean"] - .5 * (output["Profit","97.5%"] - output["Profit","2.5%"]),output["Profit","Mean"] + .5 * (output["Profit","97.5%"] - output["Profit","2.5%"]))
	output["Profit",c("2.5% NI","97.5% NI")] <- c(output["Profit","Mean NI"] - .5 * (output["Profit","97.5% NI"] - output["Profit","2.5% NI"]),output["Profit","Mean NI"] + .5 * (output["Profit","97.5% NI"] - output["Profit","2.5% NI"]))
	
	# Apply deadweight loss of taxation (from DPR AER Paper) to government spending
	output["Social Welfare",c("2.5%","97.5%")] <- output["Consumer Surplus",c("2.5%","97.5%")] + output["Profit",c("2.5%","97.5%")] - 
		1.3 * (output["Premium Subsidies",c("2.5%","97.5%")] + output["Cost Sharing Subsidies",c("2.5%","97.5%")] - 
			output["Penalties",c("2.5%","97.5%")] + output["Uncompensated Care",c("2.5%","97.5%")])
	output["Social Welfare",c("2.5% NI","97.5% NI")] <- output["Consumer Surplus",c("2.5% NI","97.5% NI")] + output["Profit",c("2.5% NI","97.5% NI")] - 
		1.3 * (output["Premium Subsidies",c("2.5% NI","97.5% NI")] + output["Cost Sharing Subsidies",c("2.5% NI","97.5% NI")] - 
			output["Penalties",c("2.5% NI","97.5% NI")] + output["Uncompensated Care",c("2.5% NI","97.5% NI")])
	
	rownames(output) <- c("Bronze P","Silver P","Gold P","Platinum P",
		"Anthem P","Blue Shield P","Health Net P","Kaiser P","Other Insurer P","HMO","PPO","Average","Average (w/ Subsidies)",
		"Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"Plan","Metal","Insurer",
		"Anthem","Blue Shield","Health Net","Kaiser","Other Insurer",
		"Bronze Cl","Silver Cl","Gold Cl","Platinum Cl","Average Cl",welfare_names)
	
	write.csv(output,"averaged_supply_output_inertia_zeroresid_variation.csv")
	
	report_shares <- FALSE
	if(report_shares) {
		for(s in 1:dim(output)[2]) {
			output[setdiff(coverage_names,"Total Coverage"),s] <- 100* output[setdiff(coverage_names,"Total Coverage"),s]/output["Total Coverage",s]
			output[setdiff(coverage_names_insurers,"Total Coverage"),s] <- 100* output[setdiff(coverage_names_insurers,"Total Coverage"),s]/output["Total Coverage",s]
		}
	}
	
	output[c("% Enrolled","Plan","Metal","Insurer"),] <- output[c("% Enrolled","Plan","Metal","Insurer"),] * 100
	
	indent_tags <- c(FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE)
	rownames(output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(output)[indent_tags],sep="")
	
	print(xtable(output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=2),
		type="latex",sanitize.text.function=function(x){x})	
	
	
# Make Counterfactual Tables for Churn Analysis and Voucher Analysis (Zero Residuals)

	basefile_columns <- c("ACA","Base","No_Inertia")
	churnfile_columns <- c("No_Churn","No_Inertia_No_Churn")
	voucherfile_columns <- c("Voucher","Voucher_No_Inertia")
	
	churn_columns <- c(basefile_columns,voucherfile_columns,churnfile_columns)
	data_churn <- as.data.frame(matrix(NA,length(fields),length(churn_columns),dimnames=list(fields,churn_columns)))
	
	for(scenario in churn_columns) {
		if(av_interaction_flag) {
			data2015 <- read.csv(paste("2015",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
			data_churn[,scenario] <- data2015
			if(scenario == "ACA") {
				data2015base <- read.csv(paste("2015",scenario,"0_supply_output_intvars_seq_avint.csv",sep=""),header=TRUE,row.names=1)
			}
		} else {
			data2015 <- read.csv(paste("2015",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2016 <- read.csv(paste("2016",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2017 <- read.csv(paste("2017",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
			data2018 <- read.csv(paste("2018",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
			data_churn[,scenario] <- (data2015 + data2016 + data2017 + data2018)/4
		}
	}
		
	
	if(av_interaction_flag & !choice_seq_flag) {
		data_churn["Profit",setdiff(churn_columns,"ACA")] <- data_churn["Profit",setdiff(churn_columns,"ACA")] -
			(data_churn["Profit","ACA"] - data2015base["Profit",1])
	}
	
	
	# Apply deadweight loss of taxation (from DPR AER Paper) to government spending
	data_churn["Social Welfare",] <- data_churn["Consumer Surplus",] + data_churn["Profit",] - 
		1.3 * (data_churn["Premium Subsidies",] + data_churn["Cost Sharing Subsidies",] - 
			data_churn["Penalties",] + data_churn["Uncompensated Care",])
	
	write.csv(data_churn,"averaged_supply_output_inertia_churn_zeroresid_intvars_avint.csv")
		
	# Extract Output
	prem_rows <- c("Bronze Premium","Silver Premium","Gold Premium","Platinum Premium",
		"Anthem Premium","Blue Shield Premium","Health Net Premium","Kaiser Premium","Other Insurer Premium",
		"HMO Premium","PPO Premium","Average Premium","Average Subsidized Premium")
	coverage_rows <- c("Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"% Switching Plan (Base)","% Switching Metal (Base)","% Switching Insurer (Base)",
		#"% Switching Plan (Setting)","% Switching Plan (Setting)","% Switching Plan (Setting)"
		"Anthem","Blue_Shield","Health_Net","Kaiser","Small_Insurer")
	#coverage_rows <- c("Bronze","Silver","Gold","Platinum","Total Coverage","Percent Enrolled","Percent Switching")
	coverage_names <- c("Bronze","Silver","Gold","Platinum")
	coverage_names_insurers <- c("Anthem","Blue Shield","Health Net","Kaiser","Other Insurer")
	claims_rows <- c("Bronze Avg Claims Resid","Silver Avg Claims Resid","Gold Avg Claims Resid","Platinum Avg Claims Resid","Average Claims Resid")
	welfare_rows <- c("Consumer Surplus","Profit","Premium Subsidies","Cost Sharing Subsidies",
				"Penalties","Uncompensated Care","Social Welfare")
	welfare_names <- c("Cons. Surplus","Profit","Prem. Subsidies","CSRs",
		"Penalties","Uncomp. Care","Social Welfare")
					
		
	output <- data_churn[c(prem_rows,coverage_rows,claims_rows),setdiff(churn_columns,"ACA")]
	rownames(output) <- c("Bronze P","Silver P","Gold P","Platinum P",
		"Anthem P","Blue Shield P","Health Net P","Kaiser P","Other Insurer P","HMO","PPO","Average","Average (w/ Subsidies)",
		"Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"Plan","Metal","Insurer",
		"Anthem","Blue Shield","Health Net","Kaiser","Other Insurer",
		"Bronze Cl","Silver Cl","Gold Cl","Platinum Cl","Average Cl")
	colnames(output) <- c("Base","(2)","(9)","(10)","(11)","(12)")
	
	
	report_shares <- TRUE
	if(report_shares) {
		for(s in 1:dim(output)[2]) {
			output[setdiff(coverage_names,"Total Coverage"),s] <- 100* output[setdiff(coverage_names,"Total Coverage"),s]/output["Total Coverage",s]
			output[setdiff(coverage_names_insurers,"Total Coverage"),s] <- 100* output[setdiff(coverage_names_insurers,"Total Coverage"),s]/output["Total Coverage",s]
		}
	}
	
	output[c("% Enrolled","Plan","Metal","Insurer"),] <- output[c("% Enrolled","Plan","Metal","Insurer"),] * 100
	
	indent_tags <- c(FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE)
	rownames(output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(output)[indent_tags],sep="")
	
	print(xtable(output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
	# Create welfare table separately
	
	welfare_output <- data_churn[welfare_rows,setdiff(churn_columns,"ACA")] * 12
	rownames(welfare_output) <- welfare_names
	colnames(welfare_output) <- c("Base","(2)","(9)","(10)","(11)","(12)")
	#welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1]) * 12
	welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1])
	indent_tags <- c(FALSE,FALSE,TRUE,TRUE,TRUE,TRUE,FALSE)
	rownames(welfare_output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(welfare_output)[indent_tags],sep="")
	
	print(xtable(welfare_output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
# Make Counterfactual tables for Robustness analysis
	
	setwd("C:/Users/eas24f/OneDrive - Florida State University/Research/Inertia RESTAT/Data")
	
	oligopoly_scenarios <- c("ACA","Base","No_Inertia")
	robust_scenarios <- c("Robust","No_Inertia_Robust")
	inertia_columns <- c(oligopoly_scenarios,robust_scenarios)
	
	example_data <- read.csv("2015Base0_supply_output.csv",header=TRUE,row.names=1)
	fields <- rownames(example_data)
	data_robust <- as.data.frame(matrix(NA,length(fields),length(inertia_columns),dimnames=list(fields,inertia_columns)))
	
	for(scenario in oligopoly_scenarios) {
		data2015 <- read.csv(paste("2015",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2016 <- read.csv(paste("2016",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2017 <- read.csv(paste("2017",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2018 <- read.csv(paste("2018",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data_robust[,scenario] <- (data2016 + data2017 + data2018)/3
		data_robust[,scenario] <- data2015
	}
	
	for(scenario in robust_scenarios) {
		data2016 <- read.csv(paste("2016",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2017 <- read.csv(paste("2017",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2018 <- read.csv(paste("2018",scenario,"0_supply_output.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data_robust[,scenario] <- (data2016 + data2017 + data2018)/3
		data_robust[,scenario] <- data2016
	}
	
	# Apply deadweight loss of taxation (from DPR AER Paper) to government spending
	data_robust["Social Welfare",] <- data_robust["Consumer Surplus",] + data_robust["Profit",] - 
		1.3 * (data_robust["Premium Subsidies",] + data_robust["Cost Sharing Subsidies",] - 
			data_robust["Penalties",] + data_robust["Uncompensated Care",])
	
	write.csv(data_robust,"averaged_supply_output_robust_zeroresid.csv")

	# Extract Output
	prem_rows <- c("Bronze Premium","Silver Premium","Gold Premium","Platinum Premium",
		"Anthem Premium","Blue Shield Premium","Health Net Premium","Kaiser Premium","Other Insurer Premium",
		"HMO Premium","PPO Premium","Average Premium","Average Subsidized Premium")
	coverage_rows <- c("Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"% Switching Plan (Base)","% Switching Metal (Base)","% Switching Insurer (Base)",
		#"% Switching Plan (Setting)","% Switching Plan (Setting)","% Switching Plan (Setting)"
		"Anthem","Blue_Shield","Health_Net","Kaiser","Small_Insurer")
	#coverage_rows <- c("Bronze","Silver","Gold","Platinum","Total Coverage","Percent Enrolled","Percent Switching")
	coverage_names <- c("Bronze","Silver","Gold","Platinum")
	coverage_names_insurers <- c("Anthem","Blue Shield","Health Net","Kaiser","Other Insurer")
	claims_rows <- c("Bronze Avg Claims Resid","Silver Avg Claims Resid","Gold Avg Claims Resid","Platinum Avg Claims Resid","Average Claims Resid")
	welfare_rows <- c("Consumer Surplus","Profit","Premium Subsidies","Cost Sharing Subsidies",
				"Penalties","Uncompensated Care","Social Welfare")
	welfare_names <- c("Cons. Surplus","Profit","Prem. Subsidies","CSRs",
		"Penalties","Uncomp. Care","Social Welfare")
					
		
	output <- data_robust[c(prem_rows,coverage_rows,claims_rows),setdiff(inertia_columns,"ACA")]
	rownames(output) <- c("Bronze P","Silver P","Gold P","Platinum P",
		"Anthem P","Blue Shield P","Health Net P","Kaiser P","Other Insurer P","HMO","PPO","Average","Average (w/ Subsidies)",
		"Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"Plan","Metal","Insurer",
		"Anthem","Blue Shield","Health Net","Kaiser","Other Insurer",
		"Bronze Cl","Silver Cl","Gold Cl","Platinum Cl","Average Cl")
	colnames(output) <- c("Base","(2)","(13)","(14)")
	
	
	report_shares <- TRUE
	if(report_shares) {
		for(s in 1:dim(output)[2]) {
			output[setdiff(coverage_names,"Total Coverage"),s] <- 100* output[setdiff(coverage_names,"Total Coverage"),s]/output["Total Coverage",s]
			output[setdiff(coverage_names_insurers,"Total Coverage"),s] <- 100* output[setdiff(coverage_names_insurers,"Total Coverage"),s]/output["Total Coverage",s]
		}
	}
	
	output[c("% Enrolled","Plan","Metal","Insurer"),] <- output[c("% Enrolled","Plan","Metal","Insurer"),] * 100
	
	indent_tags <- c(FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE)
	rownames(output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(output)[indent_tags],sep="")
	
	print(xtable(output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
	# Create welfare table separately
	
	welfare_output <- data_robust[welfare_rows,setdiff(inertia_columns,"ACA")] * 12
	rownames(welfare_output) <- welfare_names
	colnames(welfare_output) <- c("Base","(2)","(13)","(14)")
	#welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1]) * 12
	welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1])
	indent_tags <- c(FALSE,FALSE,TRUE,TRUE,TRUE,TRUE,FALSE)
	rownames(welfare_output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(welfare_output)[indent_tags],sep="")
	
	print(xtable(welfare_output,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
### Network Analysis

	# Make switching cost table
		
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=29,region = 'J98:N121',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("Base","(1)","(2)","(3)")
								
		print(xtable(data,digits=0,align='lcccc',
			caption="Sensitivty of Annual Switching Costs to Network Variables"),
			type="latex",sanitize.text.function=function(x){x})	
		
		
	# Make Counterfactual table
	
	setwd("C:/Users/eas24f/OneDrive - Florida State University/Research/Inertia RESTAT/Data")
	
	basefile_columns <- c("ACA","Base","No_Inertia")
	netfile_columns <- c("ACA_Net","No_Inertia_Exceptnet","No_Inertia_Net")
	net_columns <- c(basefile_columns,netfile_columns)
	
	example_data <- read.csv("2015Base0_supply_output.csv",header=TRUE,row.names=1)
	fields <- rownames(example_data)
	data_net <- as.data.frame(matrix(NA,length(fields),length(net_columns),dimnames=list(fields,net_columns)))
	
	for(scenario in net_columns) {
		data2015 <- read.csv(paste("2015",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2016 <- read.csv(paste("2016",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2017 <- read.csv(paste("2017",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data2018 <- read.csv(paste("2018",scenario,"0_supply_output_intvars_avint.csv",sep=""),header=TRUE,row.names=1)[,1]
		#data_net[,scenario] <- (data2015 + data2016 + data2017 + data2018)/4
		data_net[,scenario] <- data2015
	}
		
	
	# Apply deadweight loss of taxation (from DPR AER Paper) to government spending
	data_net["Social Welfare",] <- data_net["Consumer Surplus",] + data_net["Profit",] - 
		1.3 * (data_net["Premium Subsidies",] + data_net["Cost Sharing Subsidies",] - 
			data_net["Penalties",] + data_net["Uncompensated Care",])
	
	write.csv(data_net,"averaged_supply_output_inertia_net_zeroresid.csv")
		

	# Extract Output
	prem_rows <- c("Bronze Premium","Silver Premium","Gold Premium","Platinum Premium",
		"Anthem Premium","Blue Shield Premium","Health Net Premium","Kaiser Premium","Other Insurer Premium",
		"HMO Premium","PPO Premium","Average Premium","Average Subsidized Premium")
	coverage_rows <- c("Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"% Switching Plan (Base)","% Switching Metal (Base)","% Switching Insurer (Base)",
		#"% Switching Plan (Setting)","% Switching Plan (Setting)","% Switching Plan (Setting)"
		"Anthem","Blue_Shield","Health_Net","Kaiser","Small_Insurer")
	#coverage_rows <- c("Bronze","Silver","Gold","Platinum","Total Coverage","Percent Enrolled","Percent Switching")
	coverage_names <- c("Bronze","Silver","Gold","Platinum")
	coverage_names_insurers <- c("Anthem","Blue Shield","Health Net","Kaiser","Other Insurer")
	claims_rows <- c("Bronze Avg Claims Resid","Silver Avg Claims Resid","Gold Avg Claims Resid","Platinum Avg Claims Resid","Average Claims Resid")
	welfare_rows <- c("Consumer Surplus","Profit","Premium Subsidies","Cost Sharing Subsidies",
				"Penalties","Uncompensated Care","Social Welfare")
	welfare_names <- c("Cons. Surplus","Profit","Prem. Subsidies","CSRs",
		"Penalties","Uncomp. Care","Social Welfare")
					
		
	output <- data_net[c(prem_rows,coverage_rows,claims_rows),c("ACA_Net","No_Inertia_Net","No_Inertia_Exceptnet")]
	rownames(output) <- c("Bronze P","Silver P","Gold P","Platinum P",
		"Anthem P","Blue Shield P","Health Net P","Kaiser P","Other Insurer P","HMO","PPO","Average","Average (w/ Subsidies)",
		"Total Coverage","% Enrolled","Bronze","Silver","Gold","Platinum",
		"Plan","Metal","Insurer",
		"Anthem","Blue Shield","Health Net","Kaiser","Other Insurer",
		"Bronze Cl","Silver Cl","Gold Cl","Platinum Cl","Average Cl")
	colnames(output) <- c("(15)","(16)","(17)")
	
	
	report_shares <- TRUE
	if(report_shares) {
		for(s in 1:dim(output)[2]) {
			output[setdiff(coverage_names,"Total Coverage"),s] <- 100* output[setdiff(coverage_names,"Total Coverage"),s]/output["Total Coverage",s]
			output[setdiff(coverage_names_insurers,"Total Coverage"),s] <- 100* output[setdiff(coverage_names_insurers,"Total Coverage"),s]/output["Total Coverage",s]
		}
	}
	
	output[c("% Enrolled","Plan","Metal","Insurer"),] <- output[c("% Enrolled","Plan","Metal","Insurer"),] * 100
	
	indent_tags <- c(FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE,
		FALSE,FALSE,FALSE,FALSE,FALSE)
	rownames(output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(output)[indent_tags],sep="")
	
	print(xtable(output,align='lC{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
	# Create welfare table separately
	
	welfare_output <- data_net[welfare_rows,c("ACA_Net","No_Inertia_Net","No_Inertia_Exceptnet")] * 12
	rownames(welfare_output) <- welfare_names
	colnames(welfare_output) <- c("(15)","(16)","(17)")
	#welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1]) * 12
	welfare_output[welfare_names,] <- (welfare_output[welfare_names,] - welfare_output[welfare_names,1])
	indent_tags <- c(FALSE,FALSE,TRUE,TRUE,TRUE,TRUE,FALSE)
	rownames(welfare_output)[indent_tags] <- paste("\\","hskip .5cm ",rownames(welfare_output)[indent_tags],sep="")
	
	print(xtable(welfare_output,align='lC{1.7cm}C{1.7cm}C{1.7cm}',digits=0),
		type="latex",sanitize.text.function=function(x){x})	
	
	# Demand Parameter Estimates of for Network Analysis
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=4,region = 'A3:M12',header=FALSE)
		covariate_names <- data[,1]
		data <- data[,2:ncol(data)]
		
		# Create fake regression data
			
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
			
		stargazer(fake,fake,fake,fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1],data[,4],data[,7],data[,10]),
			se=list(data[,2],data[,5],data[,8],data[,11]),
			p=list(data[,3],data[,6],data[,9],data[,12]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Model Parameters",
			dep.var.caption="Estimated Model Parameters")		
			
	# Sensitivity of Demand Parameter Estimates to How We Treat Missing Values
	
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=4,region = 'A18:J27',header=FALSE)
		covariate_names <- data[,1]
		data <- data[,2:ncol(data)]
		
		# Create fake regression data
			
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
			
		stargazer(fake,fake,fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1],data[,4],data[,7]),
			se=list(data[,2],data[,5],data[,8]),
			p=list(data[,3],data[,6],data[,9]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Model Parameters",
			dep.var.caption="Estimated Model Parameters")		
					
	# Network Summary Statistics

		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=15,region = 'A1:F7',header=TRUE,rownames=1)
		colnames(data) <- c("Network Breadth","Average Network Overlap","Pairwise Overlap",
			"Pairwise Overlap, Continuing Enrollees","Pairwise Overlap, Continuing Enrollees (Excluding Previous Choice)")

		print(xtable(data,align='lC{3cm}C{3cm}C{3cm}C{3cm}C{3cm}',digits=1),
			type="latex",sanitize.text.function=function(x){x})	
	






	

# Make Regression Parameter Table (Choice Sequence)		
		
	#data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=11,region = 'A3:D20',header=FALSE)
	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=11,region = 'A4:G24',header=FALSE)
	covariate_names <- data[,1]
	data <- data[,2:ncol(data)]
	
		# Create fake regression data
		
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
		
		#stargazer(fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
		#	model.numbers=TRUE,#column.labels=c("California","Washington"),
		#	coef=list(data[,1]),
		#	se=list(data[,2]),
		#	p=list(data[,3]),
		#	dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
		#	caption="Estimated Model Parameters",
		#	dep.var.caption="Estimated Model Parameters")		
		
		stargazer(fake,fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1],data[,4]),
			se=list(data[,2],data[,5]),
			p=list(data[,3],data[,6]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Model Parameters",
			dep.var.caption="Estimated Model Parameters")		
	
# Make Regression Parameter Table (No Choice Sequence)		
	
	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=11,region = 'A55:G75',header=FALSE)
	covariate_names <- data[,1]
	data <- data[,2:ncol(data)]
	
		# Create fake regression data
		
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
		
		#stargazer(fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
		#	model.numbers=TRUE,#column.labels=c("California","Washington"),
		#	coef=list(data[,1]),
		#	se=list(data[,2]),
		#	p=list(data[,3]),
		#	dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
		#	caption="Estimated Model Parameters",
		#	dep.var.caption="Estimated Model Parameters")		
		
		stargazer(fake,fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1],data[,4]),
			se=list(data[,2],data[,5]),
			p=list(data[,3],data[,6]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Model Parameters",
			dep.var.caption="Estimated Model Parameters")		
			
	
# Make First Stage Parameter table
	
	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=30,region = 'A3:D15',header=FALSE)
	covariate_names <- data[,1]
	data <- data[,2:ncol(data)]
	
		# Create fake regression data
		
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
		
	
	stargazer(fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1]),
			se=list(data[,2]),
			p=list(data[,3]),
			dep.var.labels.include=FALSE,omit.stat=c("rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated First Stage Parameters",
			dep.var.caption="Estimated First Stage Parameters")		
			
	
# Make Demand Parameter Table

	#data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=14,region = 'A3:K42',header=FALSE)
	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=14,region = 'M3:AF65',header=FALSE)
	covariate_names <- data[,1]
	indent_tags <- as.vector(data[,ncol(data)])
	data <- data[,2:(ncol(data)-1)]
	
	covariate_names[indent_tags] <- paste("\\","hskip .5cm ",covariate_names[indent_tags],sep="")
	
	
		# Create fake regression data
		
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
		
		stargazer(fake,fake,fake,fake,fake,fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1],data[,4],data[,7],data[,10],data[,13],data[,16]),
			se=list(data[,2],data[,5],data[,8],data[,11],data[,14],data[,17]),
			p=list(data[,3],data[,6],data[,9],data[,12],data[,15],data[,18]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Demand Parameters",
			dep.var.caption="Estimated Demand Parameters")	

# Make Random Parameters Table

	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=17,region = 'A3:K42',header=FALSE)
	covariate_names <- data[,1]
	indent_tags <- as.vector(data[,ncol(data)])
	data <- data[,2:(ncol(data)-1)]
	
	covariate_names[indent_tags] <- paste("\\","hskip .5cm ",covariate_names[indent_tags],sep="")
	
	
		# Create fake regression data
		
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
		
		stargazer(fake,fake,fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1],data[,4],data[,7]),
			se=list(data[,2],data[,5],data[,8]),
			p=list(data[,3],data[,6],data[,9]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Demand Parameters",
			dep.var.caption="Estimated Demand Parameters")	
		
	
# Demographics Table

	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=1,region = 'A1:G57',header=TRUE,rownames=1)
	
	indent_tags <- !is.na(data[,2])
	indent_tags[1] <- FALSE

	rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
	data[7:(nrow(data)),] <- round(data[7:(nrow(data)),]  * 100,digits=1)
	colnames(data) <- c("2014","2015","2016","2017","2018","Overall")
	hlines <- c(-1,0,nrow(data))
	print(xtable(data,align='lcccccc',caption="Enrollee Demographics By Year",digits=0),
		hline.after = hlines,type="latex",scalebox='0.8',sanitize.text.function=function(x){x})

# New/Renewal Premiums

	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=22,region = 'A1:G6',header=TRUE,rownames=1)
	data[,c(3,6)] <- 100 * data[,c(3,6)]
	print(xtable(data,align='lcccccc',caption="New/Renewal Premiums",digits=2),
		type="latex",scalebox='0.8',sanitize.text.function=function(x){x})

# Transition Tables

	# Read in output 
		
		# Metals
		
		regression <- readWorksheetFromFile("Paper Tables.xlsx",sheet=2,region = 'A15:J24',header=TRUE,rownames=1)
	
		print(xtable(round(regression*100,1),digits=1,align='lC{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}',
			caption="Metal Transition Matrix"),
			type="latex",sanitize.text.function=function(x){x})		

		# Insurers

		regression <- readWorksheetFromFile("Paper Tables.xlsx",sheet=2,region = 'M15:V24',header=TRUE,rownames=1)
	
		print(xtable(round(regression*100,1),digits=1,align='lC{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}C{1.5cm}',
			caption="Insurer Transition Matrix"),
			type="latex",sanitize.text.function=function(x){x})		

# Elasticities and Switching Costs
		
		# Own Price Elasticities
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=29,region = 'A42:G65',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("(1)","(2)","(3)","(4)","(5)","(6)")
								
		print(xtable(data,digits=2,align='lcccccc',
			caption="Own-Price Elasticities"),
			type="latex",sanitize.text.function=function(x){x})	
		
		# Own Price Semi-Elasticities
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=29,region = 'J42:P65',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("(1)","(2)","(3)","(4)","(5)","(6)")
								
		print(xtable(data,digits=2,align='lcccccc',
			caption="Own-Price Semi-Elasticities"),
			type="latex",sanitize.text.function=function(x){x})	
		
		
		# Coverage Elasticities
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=29,region = 'A70:G93',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("(1)","(2)","(3)","(4)","(5)","(6)")
								
		print(xtable(data,digits=2,align='lcccccc',
			caption="Cross-Price Elasticities"),
			type="latex",sanitize.text.function=function(x){x})	
		
		# Coverage Semi-Elasticities
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=29,region = 'J70:P93',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("(1)","(2)","(3)","(4)","(5)","(6)")
								
		print(xtable(data,digits=2,align='lcccccc',
			caption="Cross-Price Semi-Elasticities"),
			type="latex",sanitize.text.function=function(x){x})	
		
		
		# Switching Costs
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=29,region = 'A98:G121',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("(1)","(2)","(3)","(4)","(5)","(6)")
								
		print(xtable(data,digits=0,align='lcccccc',
			caption="Switching Costs"),
			type="latex",sanitize.text.function=function(x){x})	
		
		
		
		
		
		
		# Control Function Approach
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=3,region = 'H1:M32',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("Elasticity","Semi-Elasticity","Elasticity","Semi-Elasticity","Switching Cost")
					
		print(xtable(data,align='lC{2cm}C{2cm}C{2cm}C{2cm}C{2cm}',
			caption="Elasticies and Switching Costs"),digits=1,
			type="latex",sanitize.text.function=function(x){x})	
		
		# Random Parameters
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=3,region = 'O1:T32',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("Elasticity","Semi-Elasticity","Elasticity","Semi-Elasticity","Switching Cost")
					
		print(xtable(data,align='lC{2cm}C{2cm}C{2cm}C{2cm}C{2cm}',
			caption="Elasticies and Switching Costs"),digits=1,
			type="latex",sanitize.text.function=function(x){x})	
		
		
		# Add AV and Silver Interactions
		
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=3,region = 'AD:AI32',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,5])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("Elasticity","Semi-Elasticity","Elasticity","Semi-Elasticity","Switching Cost")
					
		print(xtable(data,align='lC{2cm}C{2cm}C{2cm}C{2cm}C{2cm}',
			caption="Elasticies and Switching Costs"),digits=1,
			type="latex",sanitize.text.function=function(x){x})	
		
# Model Fit Table
			
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=31,region = 'A1:H14',header=TRUE,rownames=1)
		colnames(data) <- c("Observed","(1)","(2)","(3)","(4)","(5)","(6)")
		
		# Convert total enrollment to annual enrollment
		#years <- c(2014:2018)
		#demographic_summary = read.csv("basic_demographics_inertia_twochoice.csv",header=TRUE,row.names=1)
		#observed_enrollment <- demographic_summary["Total Enrollment","Overall"]/length(years)
		#adjustment_factor <- observed_enrollment/data["Total Enrollment","Observed"]
		#data["Total Enrollment",] <- adjustment_factor * data["Total Enrollment",]
		
		share_rows <- c("AV","Bronze","Silver","Gold","Platinum","Anthem","Blue Shield","Health Net","Kaiser","Other Insurer","Switching Rate")
		data[share_rows,] <- 100 * data[share_rows,]
		
		print(xtable(data,align='lC{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}C{1.7cm}',digits=1,
			caption="Assessing Model Fit"),
			type="latex",sanitize.text.function=function(x){x})	
		
# Demand sensitivty with network variables
		
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=4,region = 'A3:M9',header=FALSE)
		covariate_names <- data[,1]
		data <- data[,2:ncol(data)]
	
		# Create fake regression data
		
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
		
		stargazer(fake,fake,fake,fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1],data[,4],data[,7],data[,10]),
			se=list(data[,2],data[,5],data[,8],data[,11]),
			p=list(data[,3],data[,6],data[,9],data[,12]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Model Parameters",
			dep.var.caption="Estimated Model Parameters")		
			
			

# Inattention Tables
			
		# Elasticities 	
			
		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=29,region = 'AQ4:AU27',header=TRUE,rownames=1)
	
		indent_tags <- !is.na(data[,2])
		indent_tags[1] <- FALSE

		rownames(data)[indent_tags] <- paste("\\","hskip .5cm ",rownames(data)[indent_tags],sep="") 
		colnames(data) <- c("(1)","(2)","(3)","(4)")
								
		print(xtable(data,digits=2,align='lcccc',
			caption="Own-Price Elasticities"),
			type="latex",sanitize.text.function=function(x){x})	
		
		# Demand Parameters

		data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=33,region = 'A3:E58',header=FALSE)
		covariate_names <- data[,1]
		indent_tags <- as.vector(data[,ncol(data)])
		data <- data[,2:(ncol(data)-1)]
		
		covariate_names[indent_tags] <- paste("\\","hskip .5cm ",covariate_names[indent_tags],sep="")
		
		
			# Create fake regression data
			
			fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
			fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
			
			stargazer(fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
				model.numbers=TRUE,#column.labels=c("California","Washington"),
				coef=list(data[,1]),
				se=list(data[,2]),
				p=list(data[,3]),
				dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
				caption="Estimated Demand Parameters",
				dep.var.caption="Estimated Demand Parameters")	

	# Make Regression Supply Parameter Table		
		
	#data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=11,region = 'A3:D20',header=FALSE)
	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=35,region = 'A4:D17',header=FALSE)
	covariate_names <- data[,1]
	data <- data[,2:ncol(data)]
	
		# Create fake regression data
		
		fake_data <- data.frame(matrix(rnorm(100*(nrow(data))),100,nrow(data)))
		fake <- lm(fake_data[,1] ~ . - 1 ,data=fake_data)
		
		#stargazer(fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
		#	model.numbers=TRUE,#column.labels=c("California","Washington"),
		#	coef=list(data[,1]),
		#	se=list(data[,2]),
		#	p=list(data[,3]),
		#	dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
		#	caption="Estimated Model Parameters",
		#	dep.var.caption="Estimated Model Parameters")		
		
		stargazer(fake,style="aer",model.names=FALSE,covariate.labels=covariate_names,
			model.numbers=TRUE,#column.labels=c("California","Washington"),
			coef=list(data[,1]),
			se=list(data[,2]),
			p=list(data[,3]),
			dep.var.labels.include=FALSE,omit.stat=c("f","rsq","adj.rsq","n","res.dev","ser"),
			caption="Estimated Model Parameters",
			dep.var.caption="Estimated Model Parameters")		


# GCF table

	data <- readWorksheetFromFile("Paper Tables.xlsx",sheet=34,region = 'A1:H20',header=TRUE,rownames=1)
	colnames(data) <- c(2014:2018)
	
	print(xtable(data,digits=3,align='lccccc',caption="CMS Geographic Cost Factors"),
			type="latex",sanitize.text.function=function(x){x})	
		



	
	
	
	
	
	